Rippling patterns in aggregates of myxobacteria 
arise from cell-cell collisions 
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Experiments with myxobacterial aggregates reveal standing waves called rippling patterns. Here, 
these structures are modelled with a simple discrete model based on the interplay between migration 
and collisions of cells. Head-to-head collisions of cells result in cell reversals. To correctly reproduce 
the rippling patterns, a refractory phase after each cell reversal has to be assumed, during which 
further reversal is prohibited. The duration of this phase determines the wavelength and period of 
the ripple patterns as well as the reversal frequency of single cells. 

PACS numbers: 5.65.+b; 87.18. Ed; 87.18.Hf 



Pattern formation in aggregates of bacteria and amoe- 
bae is a widely observed phenomenon KM- Complex 
colonial patterns of spots, stripes and rings are pro- 
duced by E. coli 0]. Bacillus subtilis exhibits branch- 
ing patterns during colonial growth M . Two- and three- 
dimensional circular waves and spirals have been ob- 
served during aggregation of the eukaryotic slime mold 
Dictyostelium discoideum [5-8], often closely resembling 
patterns in chemical reactions P]. These patterns 
are typically modelled with continuous reaction-diffusion 
equations describing the spatiotemporal evolution of the 
cell density and concentrations of chemoattractants and 
nutrients, e. g. in the cases of Dictyostelium discoideum 
|@,|o|jl|] and E. coli ||l|,0. An alternative approach 
uses discrete models describing the motion of individual 
cells and has been previously employed to study swarm 
behavior [14-16]. 

The prokaryotic soil bacterium Myxococcus xanthus 
[0,18[ is one of the most intriguing examples for mor- 
phogenesis and pattern formation. Like the slime mold 
Dictyostelium discoideum, M. xanthus exhibits social be- 
haviour and a complex developmental cycle. As long as 
there is sufficient food supply, vegetative cells feed on 
other bacterial species, grow and divide. But when nu- 
trients run short, bacteria start to aggregate and finally 
build a multicellular structure, the fruiting body. In or- 
der to maintain this life cycle intercellular communication 
is essential. Although a phospholipid chemoattractant 
has been identified [Q recently, the key role is ascribed 
to interactions occuring through cell-cell contact during 
collisions |ll7|] . 

An experimental example for the rippling phenomenon 
|P0| , pl| is displayed in Fig. 1. Bacteria organize into 
equally spaced ridges (dark regions) that are separated by 
regions with low cell density (light regions); for a movie 
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FIG. 1. (a) Snapshot from a rippling sequence, taken from 
a time-lapse movie. White bar: 0.5 mm. (b) Space-time plot 
of the density profile along the white line in (a). 



see 1 22 . Rippling patterns were first discovered by one 
of us (H. Reichenbach) and originally named oscillatory 
waves. We examine the temporal dynamics of the density 
profile along a one-dimensional cut indicated by the white 
line in Fig. la. The resulting space-time plot is shown 
in Fig. lb and reveals a periodically oscillating standing 
wave pattern superimposed by spatiotemporal noise. 

In the following, we present a model based on the dy- 
namics of individual cells for the formation of ripple pat- 
terns during the aggregation of myxobacteria. We will 
show how certain collision rules between cells on the "mi- 
croscopic" level lead to the observed macroscopic pattern 
and reproduce the characteristics of single cell behaviour. 
The basic rules of the model are derived from experi- 
mental results by Sager and Kaiser [E3| : Within the rip- 
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FIG. 2. (a) Exemplary configuration of the model, the 
cell orientation is indicated by arrows, (b) The interaction 
neighborhood is a five nodes cross in the y, z-plane at that 
a;-position the cell is directed to (here the cell orientation is 
the +a;-direction) . 



pling phase cells are found to move on linear paths par- 
allel to each other about a distance of one wavelength. 
When two opposite moving cells collide head-on, they re- 
verse their gliding direction due to exchange of a small, 
membrane-associated protein called C-factor. Further- 
more, the model assumes a refractory phase in which 
cells can not respond to the signal. This additional in- 
gredient is necessary for the formation of the rippling 
pattern from a random configuration. The characteristic 
wavelength and period is determined by the duration r of 
the refractory phase, r is the only adjustable parameter; 
the experimental data are reproduced with a refractory 
period of five minutes. 

In the rippling phase, cells are densely packed and elon- 
gated parallel to each other. The organisation of cells into 
the aligned state can be modeled using particle-based 
models |15| , |l6| . Our model is defined on a fixed square 
grid in the x — y-plane and assumes discrete space co- 
ordinates, analagous to cellular automaton models. The 
number of cells in the model is conserved and character- 
ized by the average number n of cells per lattice site; 
cell death and replication of cells are neglected. The 
discrete z— discrete coordinate describes the number of 
cells piled up on top of each other in a given lattice 
point in the x — y-plane. The temporal update is done 
in discrete time steps. Cells mainly move along linear 
paths in the x-direction. Thus, even a two-dimensional 
model is a reasonable first step. Nevertheless, we also 
study the three-dimensional case because of its exper- 
imental relevance. The movement of individual bacte- 
ria is restricted to sheets with fixed y-coordinate. The 
coupling in the y-direction is solely due to interaction. 
A single cell is thus described by a three-dimensional 
space coordinate {x,y,z) and an orientation variable 
(f) e { — 1 (left), 1 (right)} referring to the gliding direc- 
tion. Cells interact only via head-on collisions, i. e. cells 
only sense counterpropagating cells in a certain interac- 
tion neighborhood. 

The sensitivity of a bacterium to C-factor is described 
by a clock variable ly. When a sensitive cell collides head- 
on with other cells, (the meaning of collisions will be 
specified below), it reverses its gliding direction and is 
refractory for the next r — 1 time steps, v measures the 
time since the last reversal, thus a cell with v < t is 




FIG. 3. (a) Simulation snapshot in a system of size 
220x150 with a refractory phase t — 5 (compare to Fig. 1). 
Typical pattern emerging from random initial conditions af- 
ter ca. 500 time steps, (b) Space-time plot of the density 
profile along the white line. The gray scale expresses to the 
height of the cell column on each site (black corresponds to 
high columns). 



insensitive to C-factor. Overhangs and holes in the bulk 
are prohibited. A bacterial cell is assumed to cover one 
node of the lattice, this determines the lattice constants 
of 10 /im in x- and 1 /im in y- and z— directions and the 
time constant of 1 min. 

The temporal update of the model consists of a migra- 
tion and an interaction step. In the asynchronous mi- 
gration step cells move according to their orientation to 
the neighboring site in a;-direction. If this site is already 
occupied, the cell pushes its way between cells of the 
adjacent column increasing its height by 1. With equal 
probability the cell slips beneath or above the blocking 
cell. This random process causes internal noise. There is 
also a diffusion-like noise contribution, because cells are 
assumed to rest with small probability p (in the simula- 
tions below p = 0.05 is used). Interaction takes place si- 
multaneously; every sensitive cell (v > r) checks a neigh- 
borhood of five nodes depending on its orientation (j) (Fig. 
2b). If a cell encounters at least one cell with opposite 
orientation in this neighborhood (collision), the cell re- 
verses orientation (0 -^ —(/)) and will be refractory for 
r time steps. The cell is insensitive to neighbors but 
it can still cause the reversal of other cells. Random ini- 
tial conditions and periodic boundary conditions arc used 
throughout. 

The model described in the previous section repro- 
duces the experimentally observed ripple patterns, see 
Fig. 3. The rippling pattern and the temporal evolution 
obtained in the model (Fig. 3) are in good agreement 
with the experimental data of Fig. 1. Waves propagate 



equally in both directions, their superposition forms a 
standing wave. The wavelength of the ripple pattern in- 
creases with the duration r of the refractory phase. The 
discrete model enables us to track the single cell behav- 
ior. Typically, cells move over a distance of about half a 
wavelength before they reverse their orientation (see Fig. 
4c). The pattern is easily recognizable - wavelength and 
period of the ripples have been reproduced in several in- 
depedent runs and depend only weakly on the number of 
cells in the aggregate (variations of the average number n 
of cells per lattice point between 2 and 10 do not produce 
significant changes, results presented here are for n = 3). 
Movies of the simulated rippling patterns can be found 
in ||]. 

More insight is obtained by deriving a mean-field the- 
ory of the discrete model in 2d. Such a description uses 
a hierarchy of rate equations in discrete time and space, 
which replace the discrete state variables by their aver- 
age numbers. The mean-field scheme [Q leads to the 
following set of 2r equations: 

ri{x,t + 1) = ri{x - l,t) - fr{x - l,t) + rr{x - l,t) 
r2ix,t + l)^fiix-l,t) 
r^{x,t+ 1) = r'2(x - l,t) 



rr{x, t+l) = rr-i{x - 1, t) 

h{x,t + 1) = h{x + l,t) - fi{x + l,t) + lr-l{x + l,t) 
l2ix,t+l)^ frix+l,t) 
l3{x,t+l) =l2{x+l,t) 



lrix,t+l) ^lr-lix+l,t), 



(1) 



where ri resp. h are right- resp. left moving cells that 
can reserve, while r2,...,rr resp. 12,-.-,It denote right- 
resp. left moving cells in the various stages of the refrac- 
tory phase. The functions fr and /; describe the average 
numbers of reversals of right- resp. left moving cells. 
The actual form of the reversal function is quite compli- 
cated. Since the number of particles on one site is still 
rather small, it is not sufficient to use their mean values 
in the reversal function as would be the standard ap- 
proach for rate equations describing chemical reactions. 
Instead, one has to specify the distribution of the quanti- 
ties /i, ri around their mean values and sum over all pos- 
sible states. We have used a Poissonian distribution for 
this purpose and performed a linear stability analysis of 
the homogeneous stationary state ri(x, t) — li{x, t) — ps 
and ri{x,t) = li{x,t) = pu for i = 2,....,r. This state 
describes a flat layer of cells with equal amounts of left- 
and right-moving bacteria. The actual values of ps and 
pii depend on the parameter r and n and should obey 
PR = frix, t) = fi{x, t) and n = 2{ps + (r - 1)pr). For 
r = 5 and n = 3, we find, for example, ps ~ 0.635 



and PR ~ 0.216. A detailed derivation of the mean field 
theory will be published elsewhere [Q . 

The linear stability analysis of the rate equations (||) 
with n = 3, reveals a linear instability of this fiat 
layer state against an oscillatory instability with finite 
wavenumber for r > 4 min [ p4[ . Moreover, we obtain the 
wavenumber k with the fastest growth rate for r > 4 min 
resp. weakest damping for r < 4 min and the associated 
frequency cj as a function of the refractory time r . A com- 
parison of the corresponding wavenumber A = 27r/fc and 
period T = 27r/tj with equivalent quantities extracted 
from a Fourier analysis of simulation data of the discrete 
model shows good agreement below and near the thresh- 
old r < 4 min, see Figs. 4a, b. Above the threshold 
nonlinear effects lead to a deviation of the predictions 
from linear stability analysis. It is remarkable that the 
mode with weakest damping can be observed directly for 
T < 3 min. It indicates, that the intrinsic noise of the 
discrete model drives the system out of the linear stable 
regime and excites the modes with weakest damping. 

The wavelength of the pattern in the experiment is 
about 180 pm corresponding to 18 cell lengths. The tem- 
poral period in the experiment is found to be around 10 
min. Thus, a refractory time r « 5 min in the model 
yields the correct experimental values for the wavelength 
as well as for the period. As a third quantity, we can 
measure the average reversal frequency of the individual 
cells in the simulations taking advantage of the discrete, 
particle-based nature of our model. A typical trajectory 
of an individual cell in the model is displayed in Fig. 4c. 
Most of the time the cells in the model ride with the 
ripple crest and get reflected when two crests collide. In 
other words, while the crest form is seemingly unchanged, 
most of the cells that originally constituted the crest are 
now part of a crest propagating in the other direction. 
Occasionally a cell ,, tunnels" through and continues a 
longer way with the same crest. For the refractory time 
T — b min, the reversal frequency for a single cell in the 
model is about 0.15 reversals per cell and minute in the 
three-dimensional model and 0.1 reversals per cell and 
minute in two dimensions (see Fig. 4d) which is in the 
range of the experimentally observed frequency of 0.081 
reversals per cell and minute 123]. 

While previous experiments have not provided direct 
information about the duration of a refractory phase, re- 
cent measurements of reversal rates of myxobacteria ex- 
posed to high concentration of isolated C-factor may give 
a flrst clue. Sager and Kaiser report an increase of the 
reversal rate by a factor of 3 compared to normal aggre- 
gates and an absolute reversal rate of roughly 0.3 rever- 
sals per cell per minute pS]. This suggests a refractory 
phase between 3 and 4 minutes. The reversal rates from 
the model with a refractory phase of 5 minutes would 
increase by a factor of 1.8 for the two-dimensional and 
by a factor of 1.2 for the three-dimensional model. This 
small discrepancy between model and experiment may 
indicate that the refractory time depends on the amount 
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FIG. 4. Ripple wavelength A (a) and period T (b) versus 
refractory time r for 2d simulations (solid squares), 3d simu- 
lations (open squares) and 2d mean field theory (solid line). 
(c) Single cell track with a blow-up of the marked region, d 
Reversal frequency against refractory time in 2d simulations 
(solid squares) and 3d simulations (open squares) compared 
to experiment (dotted line) and largest possible value rmax 
(solid line). 
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of C-factor and decreases at high concentrations. 

We presented a model for the formation of ripple pat- 
terns during the aggregation of myxobacteria. The re- 
versal mechanism of cells following collisions has to be 
supplemented by a refractory phase, that specifies a min- 
imum time between subsequent reversals. The duration 
of this phase determines the wavelength and the period of 
the ripple pattern. The ,, microscopic" single cell behav- 
ior agrees well with the experiments on the reversal fre- 
quency of cells. Our study strongly suggests experiments 
with single cells to verify the refractory hypothesis and to 
elucidate its biochemical basis. Moreover, myxobacterial 
rippling provides the first example of a new mechanism 
for pattern formation, namely one mediated by migration 
and direct cell-cell interaction, that may be involved in 
selforganization processes in other multicellular systems. 
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